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For a quantum state, or classical harmonic normal mode, of a system of spatial periodicity "_R," 
Bloch character is encoded in a wavevector "if." One can ask whether this state has partial Bloch 
character "fc" corresponding to a finer scale of periodicity "r." Answering this is called "unfolding." 
A theorem is proven that yields a mathematically clear prescription for unfolding, by examining 
translational properties of the state, requiring no "reference states" or basis functions with the finer 
periodicity (r, k). A question then arises, how should one assign partial Bloch character to a state of 
a finite system? A slab, finite in one direction, is used as the example. Perpendicular components 
fcz of the wavevector are not explicitly defined, but may be hidden in the state (and eigenvector \i)). 
A prescription for extracting kz is offered and tested. An idealized silicon (111) surface is used as 
the example. Slab-unfolding reveals surface-localized states and resonances which were not evident 
from dispersion curves alone. 



I. UNFOLDING IN GENERAL 

An alloy, and various other complicated systems, may 
"inherit" approximate translational order from an under- 
lying si mple r periodic system. For example, "modulated 
crystals"'i'21have a double-layer structure, each layer with 
its own "inherent" two-dimensional translational symme- 
try. Adjacent layers have translations incommensurable 
with each other. Interlayer interaction modulates both 
translations. Computations are done in a supercell, cho- 
sen to be as close a multiple as is feasible of the two differ- 
ent inherent translations. The eigenstates of the system 
remember, in some approximation, the periodicity of one, 
or the other, or both, of the unmodulated lattices. One 
should then ask, for a given supercell eigenstate how 
closely does it correspond to either of the two inherited 
families of Bloch states \k)7 "Closeness" is measured by 
weights Wi{k) between and 1, with and 1 meaning no 
inheritance or full. Dispersions can then be plotted in an 
extended Brillouin zone of a simpler system. Definitions 
of unfolding algorithms were given Ku, Berlijn, and Lee 
(KBL) et alW, and by Popescu and Zungei'*. Many ear- 
lier papers have used related versions of unfolding'^ ^ . A 
new application to alloys contains interesting views.'^The 
discussion by Popescu and Zungei'^ provides a convenient 
notation. 

A simple type of unfolding has been used in discussing 
eigenstates of multilayer^. For example, Scamarcio et 
aZ.'^ compute phonons in (Si)m(GaAs)„ multilayers, and 
relate normal modes confined in the silicon layers to vi- 
brational normal modes unfolded into a silicon crystal 
lattice. A more closely related application (but with 
weights always or 1) is unfolding optical phonons mea- 
sured in SiC polytypes into a larger Brillouin zone of a 
simpler structure of SiC. This was done by Feldman et 
al0 and Karch et al]^. 

The aim of this paper is to generalize the definition of 



unfolding. Finite systems pose a particular challenge. A 
recipe for unfolding slab calculations is proposed. 

A. Notational preliminaries 

Assume that complicated system eigenstates are avail- 
able from computation, performed with periodic bound- 
ary conditions. The "supercell" (SC) has translation 
generators Ai [i = 1,...,3). These are usually integer 
multiples of "primitive cell" (PC) translation generators 
Si {i = 1,...,3). Upper and lower case will be used 
for SC and PC properties. The integer multiple rela- 
tion is Ai = NijCij, where the 3x3 matrix N has 
integer entries. Vectors R will always be SC transla- 
tions R = '^^TTiiAi, vectors f always PC translations 
r = niSi, and when needed, the symbol x will denote 
an arbitrary location in space. 

The SC eigenstates \KJ) have the K-spa.ce transla- 
tional symmetry f{R)\KJ) = exp{iK ■ R)\KJ), for SC 
translations T by distance R. The Bloch vectors K lie in 
the supercell Brillouin Zone (SBZ) which is the unit cell 
of the SC reciprocal lattice. Its translation generators 
are Bi , where Ai ■ Bj = 2TT6ij . Reciprocal lattice transla- 
tion generators Bi of the SC and bj of the PC are related 
by the inverse of the transpose of the integer matrix Nij , 
specifically, Bi — Vectors G will always 

be SC reciprocal space translations G — '^^rriiBi, and 
vectors g — riibi will be PC reciprocal space transla- 
tions. 

The states \KJ) are not expected to have /c-space 
translational symmetry: T{r)\KJ) ^ exp{i(j))\KJ), if r 
is a translation of the PC but not of the SC. However, 
suppose the SC is nothing other than an exact A/'-fold 
repetition of the PC (where Af = detN). Then the SC 
has hidden translational symmetry, and it would be pos- 
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sible to choose its eigenstates to have the Bloch property 
4> = k ■ r, where the Bloch vector k would lie in the 
primitive cell Brillouin zone (PBZ), being related to the 
supercell Bloch vector K hy k ^ K + G, tor some SC 
reciprocal lattice vector G. In this case, unfolding would 
be an exact simplification. 

The volumes w and of the PC and SC are related 
by fl/cj = det N = A/". There are exactly TV distinct PC 
translations by multiples of di that generate the SC from 
the PC. These wiU be labeled as r^, i = 1, . . . , A/". There 
are similarly exactly Af reciprocal space translations by 
multiples of Bi that generate the PBZ from the SBZ. 
These will be labeled as Gi, i = 1, . . . , A/". These conju- 
gate sets of discrete real and reciprocal space translations 
obey the "Fourier" relations 

1 ^ - 

Y-e^G...-, ^^(^-0)^ (1) 
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These are discrete versions of the familiar "quasi- 
continuum" Fourier relations of the /c-vectors of a nor- 
mal Brillouin zone (BZ) and the discrete translations of 
a normal crystal. The crystal's translation group was 
made finite by the mechanism of Born- von Karman peri- 
odic boundary conditions. One way to understand Eq.([2| 
is as a group representation orthogonality relation. The 
discrete translations fi indeed form a finite Abelian trans- 
lation group under the closure that obtains when the su- 
percell translations are regarded as identity operators. 



B. unfolding theorem 

Suppose we have a function ^I^ ^ with Bloch symmetry 
of the superceh: g{x + R) = exp{iK ■ R)^j^{x). The 
Bloch wavevector is an arbitrary point K in the SBZ. 
This section shows that there is a unique decomposition 
of this function, ^'^j = J2g''^k+g ^^^'^ ^ functions hav- 
ing the additional Bloch symmetry of a primitive cell, 
i^^{x + r ) = cxp(zfc • r)ip^{x), where k = K + G and 
G is any of the M Gi's. The wavevectors K + G lie in 
the PBZ, but only the G = part lies in the SBZ. The 
starting function is normalized to 1, but the par- 
tial functions '4'x+g ^'^^ r\.ot. Then a candidate invari- 
ant weight for unfolding is the norm of the partial func- 
tion, Wji{G) = J (Px\ipg^(^{x)\'^. Note that the partial 
functions obtained by unfolding are not in general eigen- 
functions of any particular Hamiltonian. The function 
being unfolded, need not be an eigenfunction; it 

only needs Bloch symmetry. The translations r of the 
PC need not have a meaningful crystallographic relation 



to the translations R of the SC, except for having some 
integer commensurability relation Ai = Nij dj . 

Start by defining an operator P{K — >• K + G), that 
operates on the SC Bloch function \E'j^, and projects out 
the component that has the additional Bloch symmetry 
K + G. 

1 ^ 

P{K ^ K + G) = j^Y. r(f,)e-^(^+«>^"' (3) 
■i—i 

where T{ri)f(x) = f{x -t- r^). To prove the projective 
property, translate by a PC translation fj the function 
obtained by operating with P on a test function 

f{fj)P{K^K + G)-i'j^{r) 

= e^^+^^-^^P{K + G)*^(f) (4) 

The projected function has the appropriate PBZ Bloch 
wavevector. The proof looks a little more trivial than it 
is. One has to be careful that it is legitimate to relabel the 
J\f PC translations fi + fj as PC translations fk, because 
often fk may He outside the SC. Is it legitimate just to 
subtract the appropriate SC translation Rk to make fk lie 
in the interior of the SC? The answer is yes, for a subtle 
reason. The same difficulty occurs in proving P^ = P, or 

P(i? ^ i? + G)2 ^ E + r;-)e"^^^+'^'-("~-+*'^) 

(5) 

The simple re-labeling, that allows the double sum in 
Eq.([5]) to be Af times the single sum in Eq.(|3|, has to 
be considered with skepticism; it is not true in general, 
but only true because the P operator can operate only 
on functions with SC Bloch wavevector K. When doing 
so, any part where fi + fj lies outside the SC has the 
following effect. This translation is rewritten as fk + Rk, 
where Rk is a SC translation. The corresponding transla- 
tion operator becomes T{fk)T{Rk)- The translation by 
Rk generates the Bloch phase exp{iK ■ Rk)- This can- 
cels a corresponding phase exp{—iK ■ Rk) in Eq.(|4]), and, 
together with the identity exp(iG ■ Rk) = 1, proves the 
required theorems. 

The final step is to show that the sum of all J\f pro- 
jectors P{K K + G) is the unit operator, or equiva- 
lently, the translation T(0) by zero. This follows easily 
from Eq.Q. Therefore, any function with any SC Bloch 
symmetry K has a unique (representation-independent) 
decomposition into functions of the higher PC trans- 
lational symmetries K + Gi. This actually does not re- 
quire having an approximate higher translational sym- 
metry. The translational cells can be arbitrarily divided 
into subcells that replicate the full crystal by additional 
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translations, and Bloch functions can be separated into 
components, each of which has some higher translational 
symmetry. 



C. Unfolding formula 

The previous section suggests a formula 

WrAG) - j d^x\P{K + G)^Rj{x)\', (6) 

i.e., that the PBZ X+G-weight is the norm of the K+G- 
projected part of the SC wavefunction, 'i'gj- Using = 
P, this is equivalent to 

WjijiG) = {KJ\P{K K + G)\KJ) 



Y^{KJ\f{T,)\KJ)e 



■i{K+Gi)-rj 



3=1 



The sum rule 



(7) 



(8) 



follows from Eq.Q. This could be used on numerical 
supercell wavefunctions 'J to indicate how similar is a 

numerical iiT-state to a k = K + Gi Bloch state of the 
underlying PC. The formula uses only the eigenfunctions 
of the SC, with no reference states needed. 



D. relation to spectral function 

One of the motivations of KBiP for unfolding was to 
facilitate comparison with angle-resolved photocmission 
spectra. These provide an approximate measure of the 
"spectral function," A{k,E), rigorously definable as 



A{k,E) = {-l/Tr)trlmG{k,E + i6). 



(9) 



The Green's function matrix G„„' (fc, E) is the time-to- 
frequency transform of —i{Tct {t)c-. ,(0)). Any com- 

plete set of Bloch functions can be used to define G and 
give a spectral function A. If the set is orthonormal, in- 
variance of the trace under unitary transfomations shows 
that the spectral function is invariant. Computations for 
a disordered alloy can be done using a "supercell," an 
artificial periodic construct that enables computation. A 
sharp X-space spectrum is computed for representative 
configurations. The sharpness is an artifact, to be de- 
stroyed by averaging over an ensemble of supercells occu- 
pied by different representative local alloy configurations. 
The result, after unfolding and averaging, is a fc-space de- 
scription. Single-particle Bloch-like states have energies 



tj^j and unfolding weights Wj^j{G). These define aver- 
age values e(fc) where k = K + G. The distribution of 
weighted unfolded states determines a spectral function 
A{k, E). Green's function language allows us to interpret 
this via a complex self-energy Y,(k,E). Green's function 
theory tends to motivate perturbative understanding of 
Ti[k,E). Non-perturbative computational construction 
has only recently been implemented;'^ Computations, 
in order to sample realistically the multiple complicated 
alloy configurations, need large supercells. A single, very 
large SC computation may approach the "self-averaging" 
limit, and yield a good approximation to A{k, E) by un- 
folding. 

Consider a particular supercell with known Bloch 
eigenstates \KJ). Let there be a primitive cell with a 
complete Bloch basis \kn). Note that the states \kn) are 
not required to be eigenstates of anything except trans- 
lations, but they are required to be complete in the space 
of translational quantum number k. To any k in the PBZ 
belongs a unique K in the SBZ. (However, to each SBZ 
K there are N PBZ fc's.) Assuming the eigenstates \KJ) 
also to be complete (in the iiT-subspace, which includes 
Bloch functions oi k = K -\- G-symmetry, for any G) , we 
can write expansions 



\kn)^^ \KJ){KJ\kf 



J 



(10) 



(11) 



where cj^ is a destruction operator for the PBZ Bloch ba- 
sis function, and similarly Cj^j for the SBZ Bloch eigen- 
state. The single-particle approximation is assumed, and 
no ensemble averaging is yet performed on the SC states. 
Therefore the Green's function of the SC is 



Gjj,{K,E) - {E-e^j)-^S. 



(12) 



An alternate Green's function is defined using the PBZ 
Bloch basis, 



Gnn'{k,E) 



dte 



^{kn'\KJ){KJ\kn){E - e^^)^^(13) 



This alternate Green's function contains all the infor- 



mation of Eq.(12|, except unfolded into the larger PBZ, 



with the full details of the exact SC eigenstates obscured 
in the matrix nature and in the complex coefficients given 
by overlap matrices. Now we can use Eq.(|9| to define a 
spectral function and corresponding weight. 



A{k,E) 



T.^KAG)5{E-e^jl (14) 
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W^j{G)^Y.\^KJ\kn)\' 



(15) 



The weight Wj^j{G) corresponds to the KBL defini- 
tion, up to their assumption that the Wannier func- 
tions have the translational symmetry of the primitive 
cell unit. If the Bloch basis \kn) is complete, the weight 
of Eq.([l5| is identical to Wj^j{G) of Eq.^. The proof 
follows (tor complete, orthonormal basis sets) by insert- 
ing ig- = Y^nm \k + Gm,n){k + Gm,n\. The completeness 
relation in the SC isT-subspace requires summing over the 
M PBZ vectors k + G that map into K. This is inserted 
on either side of the translation operator in Eq.Q. Then 
simple theorems, including Eq.([2]), prove the equivalence 
to Eq.((l5l). 



II. NOTES AND POSSIBLE APPLICATIONS 

An unfolded spectrum does not contain the full in- 
formation that was in the SC spectrum. The full PC 
Green's function Gjji (K) contains full information, but 
its trace, used to find the spectral function and the un- 
folding weights, does not. The purpose of unfolding is 
to simplify in a way that provides physical insight. The 
unfolding weight defined in Eq. ([t]) has some possible ad- 
vantages. By avoiding the need for a basis of reference 
functions in a PBZ, it encourages more general appli- 
cations. For example, phonon spectra can be similarly 
unfolded. 

Among the applications that might be imagined are 
curious ones of uncertain value that could be tried. For 
example, the spectrum of a crystal with rocksalt struc- 
ture could be unfolded from the primitive rocksalt two- 
atom cell into a smaller simple cubic one-atom cell. A 
primitive fee crystal could be regarded as a supercell of 
a simple cubic structure which had alternate cells empty, 
and could be unfolded into this hypothetical primitive 
cell. Diamond structure could be considered a supercell 
which has alternating cells with two and zero atoms, and 
could thus also be unfolded to a simple cubic cell. 
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FIG. 1. Vibrational eigenfrequencies of a diatomic chain, un- 
folded onto the monatomic Brillouin zone. The spectrum 
of the monatomic chain is also shown as a thin continuous 
line. The width of the diatomic dispersion curve indicates 
the amount of weight at the unfolded wavevector. 



normal mode spectrum has two bands (acoustic, labeled 
and optic, labeled "-I-") in the BZ -7r/2a < K < 
7r/2a. Let us regard the two-atom cell of size 2a as the 
supercell (SC), and a one-atom cell of size a as the PC. 
The PBZ has -n/a < k < n/a. In the PBZ there are 
two vectors, K and K + G (where G = 7r/a), for each 
SBZ vector K. Eq.Q gives a prescription for unfolding. 
The results are shown in Figjl] 

The algebra is based on textbook procedureJ^^l^ but 
is slightly tedious. Details are in appendix [X] Bloch's 
theorem for the SC gives a 2 x 2 dynamical matrix 
D\s), which has the form 



cosf 



sin 06-'^'' 
— COS 9 



Here fig, il^ and 6 are defined by 
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III. EXAMPLE: PHONONS IN 1 D 

As a test of the unfolding formula, consider the lattice 
normal modes of a diatomic chain in one dimension. The 
Hamiltonian is 



H = T + V; T 



F 



E 



p2 



p2 

2Mi 2M2 



= — ^ [{UnA - Wn,2)^ + (Un-1,2 " Uri,l)^] (1^) 
n 

"Atoms" are alternately of type 1 and 2; they are spaced 
evenly and interact equally with neighbors on both sides. 
The separation is a and the lattice constant is 2a. The 



nl = F 



1 1 y 4,cos^{Ka) 

mI^m^j ^ M1M2 



1 1/2 



smW = 



2Fcos{Ka)/^/MiM2 



(19) 



(20) 



The eigenfrequencies are ujj. — ft^ztill, and the weights 
turn out to be 



WK±iO) = ^ [It sine] 
WK±{7:/a) = ^[l±sme]. 



(21) 
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ing in such a case? To answer this, it is appropriate 
to look at specific examples. The two examples given 
below lead to a particular way of dealing with the prob- 
lem. One wants to exhibit the standing wave nature of 
the states in the slab's interior (referred to as "bulk" 
states.) Since we care about electron states that are 
bound inside the slab, the standing waves are analogous 
to classical waves on a string with "hard wall" bound- 
ary conditions. That is, they have the character of sine- 
wave states \kn)s = {\kn) ~ \ — kn))/y^2), rather than 
cosine waves, \kn)c, as occur with open-end boundary 
conditions. There are only half as many such states as 
there are /c-vectors in the BZ, consistent with restricting 
k = (fc|| , fc2)-vectors to kz > 0. This reduces the num- 
ber of states by 2. But we need to describe Af states. 
Hard wall eigenstates are spaced twice more densely in 
kz than those selected by periodic boundary conditions. 
This observation is one guide to how unfolding must be 
done. 



Slabs extended cyclically to facilitate unfolding, (a) 
(b) is the diatomic slab with a 
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FIG. 2 

is the simple slab of Eq 
single dangling atom on one end and a corresponding Shockley 
state, (c) is the generic slab. 



The two branches K+ and K— unfold to the wavevec- 
tors k = K and k — K + ir/a. Over most of the spec- 
trum, a normal mode belongs mostly to one or the other 
/c-point, except near the SC BZ boundary K = 7r/2a, 
where sin 61 — ?> and modes belong equally to both PBZ 
fc- vectors, as shown in Fig. ([!]). If the masses evolve to 
equality (Mi = M2), the weights properly unfold the 
spectrum into the single band lu^ = (S-F/M)! sin(fca/2)|. 
It does this because the weights Wk- are 1 in the first 
SBZ and zero in the second, and opposite for Wk,+ - 



IV. SLAB STATES UNFOLDED TO BLOCK 
STATES 

DFT calculations are often done for "slabs," cells with 
a small period in x and y directions, and a long period 
in the z direction. Examples include modeling stack- 
ing faults or artificial multilayers. A common applica- 
tion is for surfaces, where the slab contains significant 
vaccuum space, with two solid/vaccuum interfaces (sur- 
faces.) Electron chemistry causes surface atoms to re- 
lax to new ground state positions. The surface alters 
the electronic structure. Most of the resulting electronic 
eigenstates are much like the eigenstates of a bulk crys- 
tal, except modified near the surface. A few states may 
be localized at surfaces. For example, dangling bonds 
are likely to be built from localized surface states. Un- 
folding can be used to illuminate the character of the 
slab states. The presence of a vacuum layer creates a 
new situation. Is there a natural definition of unfold- 



A. Monatomic linear chain 

The first example is a finite one-dimensional chain of 
M atoms, treated in orthogonal tight-binding approxima- 
tion, with a single s-orbital per atom. The s-orbital on 
the atom is denoted \£) . The single-electron states are 
IV') = ^gCt\t)^ and the eigenstates obey H\%Ij) = E\%p), 
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. (22) 

The hopping matrix element to nearest neighbor is —t\ 
farther neighbor hopping is neglected. At the surface, the 
matrix element is enhanced by a factor 7. First, let us 
solve this for the special case 7 = 1. This can be done by 
noticing that it "inherits" standing-wave states from the 
periodic solutions of the cyclic chain with 2A/'-l- 2 atoms. 
Those are the usual traveling wave Bloch eigenstates. 



|fc) 



1 



V2A^-h2 



(23) 



with eigenenergy i?(fc) ~ —2tcos{k) and k quantized as 
27rm/(27V + 2). The 27V + 2 distinct integers m lie in 
the range —Af < k < JV + I. Consider now the special 
standing-wave states 



\k)s = 



1 



y/2JV + 2 
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V2i 



(24) 
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This is both a standing sine-wave eigenstate of the cychc 
chain (for 7=1), plus, because of fc-quantization, it 
vanishes on the special atoms I = 2Af + 2 = and £ — 
Af +1. This means that all such wave-functions have two 
distinct parts separated by nodes, and that each separate 
part is an eigenfunction of the finite chain of J\f atoms 
with 7=1, after renormalizing by \/2. That is, the finite 
A/'-atom chain with 7=1 has a complete orthonormal 
set of Af standing-wave eigenstates, 

/2 

Q(fc) = {£\k{Af, 7 = 1}> = ^ ^^^ (25) 

where k is quantized as k ^ T:m/{Af +1), and where 
the integer m is non-negative, m — 1,2, ...,A/'. This 
choice satisfies Eq. ( 22 1 . The energies of these states are 
E = E{k) — — 2tcos(fc), exactly the same as for the 
"bulk" band structure. The wavevectors with m = or 
N +1 (corresponding to A; = or tt) do not work because 
the corresponding state vanishes on every atom. There is 
a close analog to vibrational standing waves with hard- 
wall boundary conditions. 

Now consider what happens when 7 is allowed to vary. 
The simplest case is if 7 is set to zero. The two end atoms 
are now decoupled, and the remaining coupled chain has 
U -2 atoms. The orbitals > for £ = 1 and £ = A/" are 
now eigenstates of energy 0, while the "bulk" states look 
just as before, except there are fewer fc's, being quantized 
in new units nm/lAf — 1), with m — 1, 2, . . . , A/" — 2. 
The atoms are re-numbered so that 2 — > 1, etc., and the 
quantized fc's are re-defined so that q vanishes at £ — 
and £ = AA- 1. 

For all other values of 7, besides and 1, careful inspec- 
tion shows that there continue to be standin g w ave eigen- 
states. However, fc-quantization is no longer^ in simple 
units of TT divided by something like Af. The actual evo- 
lution of the standing wave fc's, for a 7-atom chain, is 
shown in Fig. [s] There is a critical value, 7c — V2, 
beyond which pairs of bulk standing-wave solutions dis- 
appear, and new surface states emerge. It is simplest to 
describe if the origin is taken in the middle of the chain. 
Then the standing waves are either even (cosine waves) or 
odd (sine waves.) Similarly, the surface states are either 
even (cosh) or odd (sinh) combinations of modes that 
decay exponentially from both surfaces. Surface states 
appear both below the bulk band {E < —2t) and above 
{E > 2t). The ones below have complex fc-vector k — in 
and those above have k = ir — in. The value of k is 
related by a transcendental equation to the value of 7. 
Their energies are ±2t cosh k. The symmetric (cosh) so- 
lutions appear at 7^ > 2 while the odd (sinh) solutions 
appear at slightly larger values 7^ > 2/(1 — {2/{N — 1))). 

There is a lesson from this that applies to unfold- 
ing, namely, that in slab calculations, good standing- 
wave bulk-like solutions continue to exist. However, they 
are perturbed by the surface in three ways. First, k- 
quantization is twice as dense as with periodic bound- 




FIG. 3. fc-vectors (horizontal; between and tt) evolving 
as the surface hopping enhancement factor 7 (squared and 
plotted vertically) varies, for a 7-atom single s-orbital tight- 
binding chain. Simple periodicity of k is seen in the two cases 
7 = 1 and 7 = where surface enhancement is absent and 
the simplest bulk behavior is reproduced by slab standing 
waves. As the critical value 7 = \/2 is exceeded, two bulk 
standing waves disappear, and surface states appear with en- 
ergies above and below the bulk bands, symmetric between 
the two surfaces. Above a larger critical value (7 = \/3 for 
the 7-atom slab) two more bulk standing waves disappear, 
and anti-symmetric surface states appear. The dashed lines 
should be imagined moving perpendicular to the page in the 
imaginary k direction. These give the decay rates exp(— k^) 
and exp(— k(A/' + ! — £)) of the surface eigenstates. 



ary conditions. Instead of having TV solutions between 
k = — TT and k — n, there are Af, or A/" — 2, or A/" — 4 solu- 
tions between k = and k ~ n. Second, fc-quantization is 
no longer in easily predictable multiples of 2-k jA! . Third, 
as can be verified by direct computation, the amplitude 
of both bulk and surface states, on surface atoms, is re- 
duced by 1/7, when the surface hopping is enhanced by 
7. Similarly, in real crystals, the surface layers will dis- 
place, couplings will be altered, and bulk-like states will 
be altered, probably reduced in amplitude, near the sur- 
face. 



B. Diatomic linear chain 

There is another example, the diatomic chain of Al 
molecules, that provides further illumination. Atoms are 
paired. They have identical on-site energies (set to zero), 
weaker first-neighbor coupling —t to adjacent molecules, 
and stronger intramolecular — i' — —^t coupling within 
the molecule. The cyclic chain of 2Af + 2 molecules is il- 
lustrated in Fig|2](b). Bloch's theorem reduces this to a 
2x2 problem. Denoting the wavefunction coefficients as 
ui = Uk ex-p{ik£) I\JA! + 1 for the left atom of a molecule, 
and similarly with u v for the right atom, the Hamil- 
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tonian and eigenstates become 

( "fc ) " 7^ ( e-'^^/' ) ' ^^^^ 

with fc quantized as fc = 27rTO/(2A/'+2) and energy eigen- 
values Ek = ±|£;fe|, with \Ek\ = |i|(a2 + blf/^ and 
tan(/)fc = Qk/bk- The parameters are a^: = 7sinfc and 
&fc = 1 + 7 cos k. There are two bands symmetricaUy dis- 
tributed around E = 0. Once again, standing sine-wave 
eigenstates can be made, and can be chosen to vanish on 
two oppositely placed atoms, as shown in Figj2]^b). The 
sine-waves are eigenstates for the separated half-chains 
with the nodal atoms eliminated. Next to each nodal 
atom is an unpaired atom that remains. Each separate 
chain thus has Af molecules, plus a single weakly- coupled 
adatom on one side of the chain. The 2Af -I- 1 x 2Af + 1 
Hamiltonian matrix for the upper separated chain is 
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(28) 

There are Af positive /c's with 2A/' corresponding stand- 
ing sine-waves (the 2 being from the two orbitals per 
molecule) of the 2Af + 2-molecule cyclic chain. These 
generate states for the Af + 1/2-molecule finite chain. 
Specifically, there are states with energy — and 
Af with energy -t-|£'fe|, two for each intact molecule of 
the upper separated chain. But there are 2Af + 1 atoms 
on this chain, so one stat e rem ains to be found. It is a 
Shockley-type surface stat^^^Ell lying exactly at mid-gap, 
E = 0. Its eigenfunction, for the upper chain, has coef- 
ficients Vl = {—l/'yY and ui — 0. That is, it is zero on 
the right atom of each molecule, and decays with com- 
plex wavevector k = tt + i log 7. It alternates in sign and 
decays as exp(— ^log7), going to the right from the left- 
most atom {£ — 1). This remains a midgap eigenstate for 
ah 7 > 1. 



C. Slab unfolding 

How therefore should we "unfold" such slab calcula- 
tions? In the real 3d world of surfaces, slab states have 
fc's quantized in the x and y directions parallel to the 
slab, but perpendicular, only k^ = states are of inter- 
est. Any dispersion in fc^ derives from undesired, hope- 
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FIG. 4. Unfolded bands of an 11 atom tight-binding one-band 
model. Surface hopping is enhanced by a factor of 7 — 2. 
Symmetric and antisymmetric surface states are essentially 
degenerate. They occur at energies Q.32t above the top and 
below the bottom of the bulk bands at ±2t. The bulk band 
dispersion E{k) = —2tcos{k) is plotted, and agrees with the 
wavevectors found for bulk states in the slab's interior. 

fully small, interaction between the slab and its periodic 
images. These images were included only for computa- 
tional convenience, to permit the use of periodic codes. 
The true slab states are either surface states, or else 
they are extended, and hence, standing waves. Even 
though they have quantum number kz—0 from the band- 
structure code, nevertheless, the extended states relate 
to superpositions of bulk Bloch states with equal contri- 
butions from kz and —kz for some value of kz- The va- 
lence and low-lying conduction standing-wave states are 
bound in the potential well of the solid. Thus they have 
sine- wave rather than cosine- wave character. We wish to 
know the approximate relevant kz = \kz\. Their weights 
will be the same at +kz and —kz- The maximum value 
of kz is certainly tt/c where c is the (atomic scale) layer 
spacing within the slab. But the kz values are definitely 
not related to multiples of 211 jAf where Af is the number 
of layers. Instead, their average spacing is 27r/(2A/'+ 2). 
Precise spacing is only obeyed for very special model sit- 
uations. A standing wave of a thick slab should have a 
sharply defined, but slightly unpredictable, kz character- 
izing its behavior in the interior, and more complicated 
wavefunctions near the boundaries, predictable only by 
detailed calculation. 

The procedure that seems most general and sensible is 
illustrated by Figj2]jc). The slab inherits an approximate 
periodicity c in the perpendicular direction, and has Af 
well-defined layers plus possible adatoms or reconstruc- 
tion at the surface. This gives an approximate thickness 
Afc + 8. Imagine instead a cycle of length (2A/' -I- 2)c. 
If the surface is thick with adatoms, [2Af + 4)c might 
be preferable. Reflect the eigenstates of the actual slab 
through the separating midline of Fig. [2j This generates 
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a candidate continuation onto the (2A^+2)-loop. But be- 
cause of the hard-wah standing wave nature of the bulk 
states, and sine-hke disappearance of the bulk states as 
they merge into the vacuum, the reflected wavefunction 
should have its sign changed to give it sine-like charac- 
ter. Finally, use the unfolding algorithm explained above 
in Sec. II CI The basis fc-vectors can be chosen to have 
= 7rm/(A/' -I- 1) for m = 1, 2, . . . , A/". These will not 
exactly conform to the actual standing-wave eigenvector 
wavelengths, but provide a convenient basis set for un- 
folding. Surface states will unfold to a broad range of 
fcz's, while bulk-like standing waves, for thick slabs, will 
have weights distributed more narrowly. 

Finally, it is sensible for plotting E versus to use 
a horizontal bar to indicate the width in k^ that a 
state unfolds to. The mean and rms fc-values can ob- 
viously be defined as A: = Sfc ^(^)^/ X^fe ^(^)' 
= J2k W{k){k - kfl Y.k ^ik)- 

As an example. Fig. |4] shows the results for the single- 
orbital Id tight-binding chain of 11 atoms thickness. The 
surface enhancement factor was chosen as 7 = 2, large 
enough so that both even and odd surface states occur, 
both above and below the bulk bands. The layer thick- 
ness of 11 is sufficient that very little splitting occurs 
between even and odd. The figure shows by circles the 
"exact" values of k^ that characterize the interior. These 
are found by numerical solution of transcendental equa- 



tions derived from Eq. ( 22 ) . Their energies agree with the 
eigen-energy — 2tcos(fc)~^ 

V. SILICON (111) SLAB UNFOLDED 

We performed an A/" = 18 double-layer silicon slab cal- 
culation using the SIESTA^IS implementation of density- 
functional theory (DFT), and the PBE-GGAEs! exchange- 
correlation potential. The slab orientation is perpendic- 
ular to (111). Charge self-consistency was obtained using 
128 k-points (all perpendicular to the (111) axis.) A cu- 
bic lattice constant a — 5.499A (compared to experimen- 
tal a = 5. 43 A) minimizes the total energy of the 2-atom 
bulk primitive cell (PC), using the standard double-zeta- 
polarized (DZP) basis orbital set with 13 orbitals per 
silicon atom. The unit cell has 36 atoms, one per layer 
(two per double layer.) Translation vectors within the 
(111) plane are a = f(-lOl) and f(0 - 11). Since 
the aim is to study unfolding, the simplest possible sur- 
face is used, namely a completely unrelaxed cut between 
adjacent double layers. The bottom of the slab is sym- 
metrical with the top. The top-most and bottom-most 
atoms have one dangling bond each. The shortest trans- 
lation vector in the (111) direction is c = a(l 1 1). The 
vectors (a, b, cj define an hexagonal cell with 6 atoms ar- 
ranged in three double layers, as shown in FigjS] The 18 
double layers define a section of bulk with vertical dis- 
tance 6c. On top of the slab is a large vacuum layer, and 
periodic repeats with period 50c. 

The computed band structure is shown in Figj6j Note 




FIG. 5. Section through the silicon slab, seen along a — b. The 
vector c || (111) is vertical, with vacuum at the top. Vector a+ 
b, defined in the text, is horizontal. Two types of translational 
cells are shown. The parallelogram outlines a diagonal slice 
through a primitive diamond cell with two silicon atoms. The 
rectangle outlines a slice through a hexagonal cell containing 
6 atoms, arranged in three bilayers. The bilayers stack in the 
fee ABC pattern. Also shown are contours of the real part 
of ■(/i for one of the pa surface bands in the gap, shown by an 
arrow in Fig[7] Contours are evenly spaced; pale and dark 
grey indicate opposite signs. 



that the Fermi level lies in the middle of a doublet of 
surface states. The surface is metallic. It is certainly 
possible that including spin-polarized options in the DFT 
calculation would have split up and down spins, perhaps 
creating a magnetic insulator instead of a metal. It did 
not seem important to test this, because the unrelaxed 
ideal surface can only be regarded as an over-simplified 
model. 

The primitive unit cell has two atoms. A conven- 
tional description would use vectors d,b,d, where d — 
c/3 — (a -f b)/3 is shown in Fig s] as the slanted side of 
the parallelogram. It has vertical height |c/3|. Two com- 
ments can be made. First, a 2-atom hexagonal cell de- 
fined by a, b, c/3 is an allowed unit cell; it translates by 
multiples of a, b, d to fill space, like bricks, with a hor- 
izontal off-set. Second, the primitive cells fill the slab 
unit cell (supercell) by symmetry translations that slant 
relative to the c-axis. This forces simple modifications in 
applying the unfolding algorithm. 

SIESTA gives expansion coefficients, denoted here by 
{N\J), for the slab eigenstate with band index J, in terms 
of the slab local basis orbitals |A^). We examine only 
eigenstates with K = 0, and ask how their approximate 
translational symmetry unfolds onto wavevectors parallel 
to c. The slab comes from a section of bulk of thickness 
18c/3 where c = |c| = V^a. The longest standing-wave 
wavelength is then approximately Amax = 2 • 19 • (c/3). 
The factor 2 is because the longest sine-wave has only a 
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are given by 




4jt/3a (01-1) 



k-vector 



FIG. 6. Band dispersion for the (111) slab of 18 double layers, 
with 36 atoms per cell. Two directions perpendicular to (111) 
are shown. The Fermi level (dashed line) intersects a doubly- 
degenerate surface state. This state is a result of one dangling 
bond per surface atom (and two surfaces), obtained because 
atoms are frozen in perfect terminated bulk positions. The 
arrow points to occupied surface states just above the top of 
the valence band, one of which is shown in contour in Fig[5] 



half wave-length in the slab. The factor 19 is used instead 
of 18, following Sec IV The extra distance c/3 estimates 
the additional distance for a typical sine-wave state to 
vanish at the edges. The fcc-character of states should 
be approximately quantized in units k = m ■ (27r/A,iiax)j 
m going from 1 to = 18, as required to represent the 
states of the 18 double layers. 



The unfolding algorithm of Sec IV 
plemented as follows 



id KBiP is 



First, a mirror copy of the slab is 
positioned adjacent, with the extra separation c/3. Sec- 
ond, wavefunctions | J) must be extended from the orig- 
inal slab to its mirror. Third, this 2JV + 2 double-layer 
slab is repeated periodically (vacuum is neglected except 
for the two extra c/3 slices.) The mirror image slab or- 
bital N' corresponding to orbital N of the real slab is as- 
signed the expansion coefficient {N'\J) = —{N\J). The 
extra space c/3 contains a missing double-layer, whose 
2 missing atoms contain contributions to state \J) with 
{N\J) = 0. All coefficients are divided by 1/2 to main- 
tain normalization. Fourth, a mapping is constructed be- 
tween orbitals of the slab and the 26 "reference orbitals" 
that form the SIESTA basis states of the primitive 2- 
atom cell. The transverse components of the positions 
are ignored, as only the relative position in the (111) 
direction is relevant for kz unfolding. 

The slab eigenstates are unfolded onto the 2J\f + 2 
wavevectors k = (0, 0, kz) = {rmr/ {JV + l)a)(lll), where 
m = —A/", . . . , —1, 0,1, ... ,J\f + 1. Only values m = 
1, . . . ,M are needed, as confirmed below. The weights 



Wj{kz) = ^{kzn\J){J\kzn'){kzn'\kzn). 



(29) 



This is the altered version of Eq.( 15 ), the difference aris- 



ing from the non-orthogonal nature of SIESTA basis 
functions. The elements {kzn\J) are given by 



{kzn\J) 



N 



(30) 

The sum over N runs over all orbitals in the slab plus 
mirror slab construction. r(N) and 5nn(N) comprise the 
two components of the orbital mapping. n{N) is the PC 
orbital corresponding to the slab orbital N and r{N) is 
the position vector connecting n{N) to N . r\){n) is the 
position of the orbital n{N) within the PC. In the formula 
for the weight Wj{kz), {kzn'\kzn) is the overlap matrix 
element between the two Bloch basis states associated 
with the local orbitals n and n' . This overlap matrix is 
extracted from a SIESTA calculation for the bulk silicon 
PC. 

The weights for kz and —kz are found to be identi- 
cal, as expected. The weights for kz = and kz = n 
are zero, so there are 18 independent weights for each 
state I J), as expected. On average, the sum of the un- 
folding weights Wj{kz) satisfies the sum rule given 
by Eq. (8). More specifically, the average value of 
Efc, Wj{kz) = 0.998 ± 0.114. The sum rule is not ex- 
act, because of non-orthogonality of the SIESTA basis 
orbitals. The plot of E versus kz in FigjT] shows for each 
eigenstatc near the band gap, the mean and rms variance 
of the distribution of kz that the state unfolds to. 

Notice that there are three states visible in the gap in 
Figj?) each of which is actually a doublet. The position 
and width of their fc-distribution reveals that they are 
paired surface states, lying on top and bottom surfaces, 
almost completely de-coupled. The higher energy dou- 
blet involves (lll)-oriented p-states, which bond strongly 
to each other when the separated surfaces are allowed to 
recombine. Upon coupling, one state moves down to the 
top of the valence band, while the other moves up, above 
the anti-bonding s state, in the conduction band. The 
lower energy quartet, when viewed without unfolding, in 
Fig|6j would not be easily recognized as localized, except 
that the unfolded version in FigjT] makes its nature clear. 
We examined contour plots of all four of these states, to 
confirm the localized nature. One of these states is shown 
in FigjS] They are p-states oriented parallel to the plane, 
bonded weakly to neighboring surface atoms. It is not ob- 
vious that opening the vacuum between bi-layers should 
so much affect their localization. 

Of the higher-lying anti-bonding states, inside the con- 
duction band, shown in Figj?) three states near —0.6 eV 
seem nearly degenerate, and two of them show anoma- 
lous fcz-broadening. These two are actually a quartet of p 
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regarded as a piece of a periodic medium with Af + 1 
(rather than the nominal Af) vertical repetitions of the 
bulk structure, the last repetition being vacant of atoms. 
To accomodate the half-size quantization of stand- 
ing waves, the slab is then mirrored to a 27V + 2-layered 
slab with periodic boundary conditions. Wavefunctions 
on the mirror image are negatives of their mirror coun- 
terparts. Finally, these modified states are unfolded in 
the normal way, revealing the hidden character of the 
states. Localized surface states, and partially localized 
surface resonances, are revealed. Not all of these prop- 
erties were expected or discovered from the conventional 
dispersion plot (Figjoj), but were uncovered by the un- 
folding algorithm. 



FIG. 7. The k = (0, 0, 0) bands of the 18 bi-layer Si (111) slab 
(the states shown at F in Fig|6| are here unfolded to reveal 
their relation to bulk states with k = (n /a){(^, 0- Smooth 
curves show the bulk dispersion in this direction. Horizontal 
bars indicate the position and width of the peak of the distri- 
bution in ( for each eigenstate of the slab. The dashed line 
indicates the Fermi energy, and the arrow locates the surface 
state shown in Fig. [5] 



states, localized (according to contour plots) on top and 
bottom surfaces. They are anti-bonding counterparts to 
the in-gap surface state quartet discussed above. The 
third state with the smaller fc^-distribution is quite dif- 
ferent. It is a single delocalized state of s-character, a 
member of the rapidly upward-dispersing band of anti- 
bonding s-states. The localized quartet cannot be true 
surface states, since they do not occur within an energy 
window where bulk states with k ~ (0, 0, fc^) are missing. 
They are degenerate with the rapidly dispersing band, 
and are therefore surface resonances. Only a larger cell 
(more than 18 double-layers) could reveal the resonant 
aspect. Of course, none of these appear in actual exper- 
imental silicon (111) surfaces. There is an easy route to 
a "2 X 1" reconstruct ior(i^, and a more difficult route to 
the ground-state 7x7 reconstructioiP^. 



VI. SUMMARY AND CONCLUSIONS 

This paper has done two things. The notion of unfold- 
ing has been given a more general interpretation. A can- 
didate general definition is made that uses only the states 
under direct consideration, needing in principle no basis 
of Bloch states of the simpler primitive cell. Only the 
properties of the complicated states under partial trans- 
lation are used. 

The idea of unfolding is then generalized to systems 
that are finite in one direction. Specifically, we consider 
a slab, which "remembers" or "inherits" the intrinsic pe- 
riodicity c that would have been exact in bulk. The vac- 
uum part of the slab is ignored, and the occupied part 



Appendix A: 1-d phonons 

The Hamiltonian of the diatomic chain of vibrating 
atoms is Eq.(16). According to Eq.Q, the weights that 
we want for this system are 



WkAB) 



l + {K,±\T{a)\K,±)e 



-i(K+B)a 



(Al) 



The chain has periodicity 2a. The translation T{a) 
gives the "inherited" approximate translational symme- 
try. The factor exp{~iBa) is either -1-1 (when K is un- 
folded onto itself by i? = 0) or — 1 (when K is unfolded 
onto K + B = K + ir/a in the "second" BZ.) 

First, make the standard SBZ Bloch- wave (Fourier) 
transform for displacements, Un^i — >■ Ui{K). Looking for 
normal modes that oscillate sinusoidally with frequency 
ujk, Newton's laws become 



w^Mlu) = f\u) 



ui{K) 

U2{K) 



(A2) 



The mass matrix M is diagonal with elements Mij = 
It is convenient to make the mass-weighting 



transform |s) = M^/^|m), which converts Eq.(A2) to the 
standard Hermitean problem 



K\ 



D\s) 



(A3) 



where the mass-weighted force matrix D is 

/ ^ 2cos(Ka) ^-iKa \ 

D^F\ . (A4) 



2cos{Ka) ^iKa 
VA/1M2 



2 



It is convenient to rewrite this as done in Sec. |III| Eqs.(|T7[ 
20 1 The eigenvectors of the dynamical matrix (EqjlT]) are 



sin(|)e-^-^'^/2 
cos(|)e+*-^''/2 



cos(f)e-^^°/2 
— sin( 



B\ +iKa/2 



(A5) 

We now need to translate an eigenvector. To find 
the form of T{a)\si), it is helpful to display the full 
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coordinate-space vector |s) of raass- weighted displace- 
ments for a Bloch wave with wavevector K , 



( 



\ 



526 
526 



2iK{n-l)a 
2iK{n-l)a 
2iKna 
2iKna 
2i(n+l)Ka 
2iK{n+l)a 



(A6) 



V: / 

This tehs us that the rule for translation is 



T{a) 



S2{K) 



S2{K) 

51 (i^) 6+2*-^" 



(A7) 



The diagonal matrix elements of the translation by a are 
then 



(5± |r(a)|5±) = Tsinf 



(A8) 



We can now evaluate the unfolding weights, Eq.(|A4|). 
Note that in doing so by use of Eq.(A8), we make a 



somewhat arbitary choice, that the eigenvectors needed 
m Eq.Q are the mass- weighted displacements |5), rather 
than the simple displacements \u). This has an advantage 
of simplicity, but more than that, the |5) states are or- 
thonormal by virtue of the ordinary Hermitean nature of 
Eq.(A3|. The states \u) have instead, the generalized or- 

|u,/3) = 5a,i3- This would 
Rather than 



thonormality relations, (u, aj 
have altered the sum rule on the weights, 
adding to 1, they would add to a mass and if-dependent 
number. The answers for the simpler |5) eigenvectors is 
given in Eq.(21) and Fig. [ij 
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